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ABSTRACT 

We develop a correspondence between arbitrary tensors and matrices based on the use of 
Kronecker products and associated identities. Utilizing the rules of matrix differentiation we derive the 
vacuum Einstein field equations as a differential-matrix equation. This formulation may facilitate their 
efficient use in numerical relativistic models. 



The general covariant field equations of general relativity find a natural representation using 
tensors. They consist of ten independent nonlinear second-order differential equations for the metric 
tensor components g^ v in terms of the four space-time coordinates, and in a vacuum reduce to the 
vanishing of the Ricci tensor 1^ = 0, which is obtained from the only meaningful contraction of the 
Riemann curvature tensor R^va. In the last several decades many researchers have worked with the 
field equations recast using the tools of differential forms^ (the exterior calculus) and in the context of 
the null-tetrad, or Newman-Penrose, formalism. 3 Here we shall demonstrate that the vacuum field 
equations can be expressed as a system of first-order differential-matrix equations. 

Tensors are convenient tools based on indicial structures and the summation convention, and as 
such enable N-dimensional arrays to be easily manipulated via transformation between upper 
(contravariant) and lower (covariant) indices and via contraction of indices. A matrix formulation 
regards all objects as two-dimensional structures which obey strict manipulative rules. However, 
equations posed in matrix form have a substantial repertoire of estimation and analyses tools available 
which can be utilized in numerical solutions; this may be especially important in numerical relativistic 
studies. As emphasized by E.T. Bell, matrices are "...a trivial notational device (which) may be the germ 
of a vast theory having innumerable applications. "4 Matrices may be thought of as very powerful 
organizational tools blessed with wide-ranging formal mathematical insights. 

The connection between a 4kx4m matrix A, for four space-time dimensions, and tensor 
A PlPz ' Pm aia2 ... C q c , where for each ocj or (3j the subscripted index j e {1,2,3,4}, is established via the Kronecker 
product, "®", defined between matrices and vectors. 5 The formalism is easily generalized to an arbitrary 
number of dimensions, j, so A is a jkxjm matrix. For any two matrices A and B, the Kronecker product 
A ® B has every element of A multiplied by the entire matrix B, i.e., 



for an nxm matrix A. This operation greatly increases the dimension of A®B; in general the 
Kronecker product is not commutative so A ® B # B ® A. Consider a set of arbitrary four-element 
column vectors {a a i,a a2 ,- ■ -^a^bp^bp!,- ■ -,bp m } , where the subscripts designate different vectors 
corresponding to the index and not components, the dimensions of A and ordering of elements in A is 
specified by the Kronecker product ordering of these vectors in the scalar i.e., 



A ® B = 



anB a 12 B a 13 B 
a 2iB a 2 2B a23B 
a3iB a32B a33B 



ai m B 
a2mB 
a3 m B 
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We call this the \|/-representation of a matrix A . 

Interchanging any two indices in the tensor A PlPa '"'Hxioe. . at is equivalent to interchanging the 
corresponding indexed vectors in \\r and reordering the elements in A such that \|/ retains its same 
numerical value, i.e., remains invariant. Raising or lower indices in A PlP2 ^ m ai a 2 ...at corresponds to shifting 

the corresponding indexed four-vector to the proper location within the Kronecker products on each side 
of the matrix A in \\i and again the elements of A are reordered such that \\f remains invariant. The row 

vector operating on the left of A specifies the covariant indices and the column vector operating on the 
right of A specifies the contravariant indices in our convention; a covariant vector is a column vector 
and a contravariant vector is a row vector. Contraction of an upper and lower tensor index, which is a 
generalized trace operation,^ corresponds to setting equal the corresponding vectors in \|/ and taking the 
sum of appropriate elements in A . 

Interchanging two indices in a 13 ' 132 Pm ai a 2 a k and keeping \|/ invariant serves to remap either rows 
or columns of A in analogy with the elementary row and column operations. ^ For example, 
interchanging the oci and ocj covariant indices of the tensor is equivalent to the following correspondences 

(a* ® a£ 2 ® . . .® a T ® . . . ® aj. ® . . . ® a T k ) A (b Pl ® b fc ® ■ ■ ■ ® b p J 
= ® a^ 2 ® . . .® a T ® . . . ® a T ® . . . ® a T k ) A v 'j (b pi ® bp 2 ® ■ ■ ■ ® bpJ 

s, A - APlPz- ■ -Pm 

^ A — A aia2- ■ -ai - aj- ■ -ak 

v A Vij _ A PlP2- ■ -Pm 

' ' — aiaz ■ -Op Cii- ■ -ak > 

where the notation A VlJ designates a vertical row reordering interchanging those rows of A corresponding 
to those operated upon by the same indexed vectors on the left. Interchanging vectors on the right of A 
in \|/, i.e., Pi and Pj, corresponds to a horizontal reordering of columns in A and this is designated as A hlj . 
In each case the reordering of rows and/or columns in A is performed such that \|/ remains invariant. 

Furthermore, the vertical and horizontal transpose operations are commutative, and for two matrices X 
and Y and a row operation § , (|)(XY) = (|)(X)Y, and for a column operation |i, p.(XY) = X|i(Y).7 

The above correspondence relations furnishes a rigorous connection between indicial tensor 
quantities and formal matrix mathematics. When considering some general tensor A 13 '* 32 ' ^ m aia2 . . a k , from 
the point of view of doing numerical analyses, forming a array may be a straightforward matter of 
rearranging all of the components of the tensor (which may be treated as a km-dimensional matrix or 
array) into a two-dimensional field and then having the computer keep track of their disposition 
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throughout the calculation. Our approach is significantly different in that we are interested in 
formulating problems, which hitherto have been treated entirely using tensor techniques and have 
appeared to be natural to tensor mathematics, fundamentally in formal matrix structures. In other 
words, we desire to cast the mathematical form of the physical problem as a matrix equation wherein the 
unknown quantity is a matrix; such equations may be nonlinear, algebraic, or differential. An important 
property of these equations is that the matrices comprising the factors in the terms do not commute so 
that the solution of these equations are extremely challenging. Matrix equations arise in many 
disciplines, including control theory, queuing theory, 8 and stability,^ and even optics. 10 

The development of the field equations begins with the scalar quadratic form ds 2 = du T Gdu, 
where G is the symmetric 4x4 metric tensor matrix with elements g^ v - Here du T = [dm du2 dir? diu] are 
the coordinate differentials and ds 2 the line-element. Performing the standard variation along ds yields 
the geodesic equation for the coordinate vector u 

d 2 ^ + r (du ® du) = r = G -l 

ds 2 ^ds ds/ 



[dGdGdG dG] 

dui du? du^ du4 



UdGf 
2Vdu/ 



(2) 



where the vector derivative of G is defined such that the first column of dG/du is the derivative of the 
row-major ordered elements of G with respect to m, etc.H The matrix r is a 4x16 matrix analog to 
the Christoffel symbol, each 4x4 block dG/du; is the scalar derivative of G with respect to ui. In the 
derivation of (2) the chain rule formula for derivatives of matrices with respect to a vector is used^ 

d(F , GH) = (FG ® Ijdtt + (F ® H T )<^ + (I r ® H T G T )^ , (3) 
du du du du 

where F is rxs, G is sxt and H is txu, and Ij are jxj identity matrices. Here G is treated as a unique 
tensor independent of the above rules to enable the second derivative term to be isolated by allowing 
G" 1 to exist. 

The Einstein field equations are derived by taking the difference of the twice covariant 
differentiated scalar invariant / = a T (du/ds) with its commuted derivatives for some vector a and 
contracting that tensor with respect to its one contravariant and a covariant index. In our formalism the 
Riemann curvature matrix is defined as 

ds2 dsi dsi ds2 \ ds dsi ds2/ 
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where the subscripts on the ds and du are used only to preserve ordering during the computation. All 
terms consisting of the vector derivatives of the vector a cancel out. The Riemann curvature matrix R , 
in direct analogy with the tensor form R £ live, IS 



/ T\ BTv 23 / T\ BT 

r = (r T ® i G ) v V - (r T ® i G )r T + (^-j - (<^J , (4b) 

where the vertical row reordering operation "V23" is defined by the relation 

[ML <g> @ du|] T @ = (d^T duT du|\ x Ig) V23 

Ids ds 2 dsi/ \ ds dsi ds 2 / 

from which it is a straightforward matter to obtain (T T ® Ig) V2 \ and similarly the third term in (4b). For 
r = [Ti T 2 T 3 r 4 ] , which consists of 4><4 blocks H, r BT = [rj rl rl rl] ; we call T B the block-transpose of r 
whereby the individual H are not transposed; T T is the standard matrix transpose and Ig is a 4 X 4 identity 
matrix, having the same dimensions as the matrix G. Important formulae and identities used in this 
computation include Equation (3), and vec(AXB) = (B T ® A)vec(X) , where vec(X) is the column-major 
ordered form of any matrix X ,4 and the Kronecker relationships a T X(M ® b) = b T X BT (M ® a) and 
(Ig ® a T )X B = c(X b )(Ig ® a) , for conformal matrices X and M , and column vectors a and b , and c(X B ) is 
defined so that row i is [vec(Xi)] T , for a square block component Xi of X , where the blocks are arranged 
vertically within X. Note further that in Equation (2) that c[(dG/du)] is the first term in F. 

The matrix R has 64 rows and four columns corresponding to the three covariant indices and one 
contravariant index in R E ^o, as can be seen from comparing (4a) and (1). Contraction to obtain the 
Ricci tensor from R%v entails — for this choice of indices — to taking the trace of each of the sixteen 
4 X 4 blocks in R; these are the sixteen components of the vacuum field equations, and is ordered as a 
1 6 X 1 matrix, which is designated as 0. = are the vacuum field equations, and is a covariant- 
indexed quantity as it is a 16-element column vector, i.e., its \|/ -representation via (1) has vectors only 
on the left. 

In applying a matrix formulation to the determination of the Riemann curvature tensor and 
hence leading to the Einstein vacuum field equations, one may ask why we even bother to discuss tensors 
if our intent is to remain fully within one domain. The reason that the \|/-representation arises, and quite 

naturally, is because of how the field equations are developed from derivatives of a scalar invariant 
associated with parallel transport of a vector along different geodesic trajectories in a curved space 
embodied in (4a). For / = a T (du/ds) the following expressions are readily derived 

dl = du T (. r BT (lG ® a ) + 4a W 
ds ds I dulds 
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Each expression is equally valid and signifies a second-rank tensor of the mixed form A^ v in the first case 
and the covariant tensor A^ v in the second case. The transpose of each expression is also a scalar and 
corresponds to the mixed tensor A^ y , and the contravariant tensor A^ v , respectively. Interchanging the 
imaginary indexed order of the coordinate derivatives outside the brackets is in keeping with the above 
rules of index switching and hence the two indices in this second rank tensor. It is important that only 
one representation be adopted throughout a derivation, such as that leading to (4b), and the proper 
implementation of a contraction yields the system of equations which interest us. Contraction on the 
\|/-representation is a generalized trace operation with respect to the chosen vector on the left and on 

the right of A, whereby the contracted matrix has l/4th the number of rows and columns. That is, 

contraction may be performed on a rectangular matrix such as the Riemann curvature matrix R as it 
appears in (4a). For the first equation above, corresponding to A^ v , contraction is simply the trace of 

the bracketed quantity. All of this goes to show the motivations and power and flexibility of tensor 
mechanics, though the \|/-representation appears to furnish a bridge between formal matrix mathematics 

and tensor mechanics, which to our knowledge has not hitherto been clearly delineated in the physics 
literature. 

Unfortunately the highly organized and rigid rules of matrix mechanics has a number of 
drawbacks. The first is that the derivation leading to (4a) makes use of non-intuitive derivative 
expressions such as the chain rule in (3). The second is that in order to keep all terms consistent in 
regards to the \|/-representation it is necessary to use complex Kronecker matrix identities and to adopt 

usage of a number of unconventional matrix representations, i.e., the vertical and horizontal reordering 
operations. The third is that many of the familiar identities with the Riemann tensor and its derivatives 
(the Bianchi identities) as well as the definition of the Weyl tensor are not obvious or natural to the 
matrix formalism. The major advantage resides in the compactness of the resulting equations for 
numerical work; possibly new inferences about the structure of the field equations can be made from more 
detailed examination of our expressions utilizing the substantial theory developed in linear and 
multilinear algebra. One such, perhaps already quite obvious, inference is that the definition of r and the 
contracted Riemann expression form a system of matrix equations. The first is a linear differential 
equation in G and the latter a Riccati-type equation in F. Substantial work on the existence and 
determination of solutions of matrix Riccati equations exists, 13 an d as such our representation of the 
field equations may enable new mathematical techniques to be brought to bear on the struggle to crack 
the field equations for various physical problems. The matrix expression are also much easier to handle 
with regard to implementing numerical algorithms using many standard software packages (IDL, Matlab, 
and Mathematica, for example) even with our unconventional reordering operations. Conceivably the 
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matrix equation formulation may be useful for developing algorithms to tackle the equivalence problem 
in general relativity involving finding practical algorithms or analytic tests to show when two different 
metrics are identical even though they are expressed different coordinate systems. 14 

The author wishes to thank Eric Le Bigot for his very helpful discussions on this topic. 
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